{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Grid states under uncertainty\n",
    "\n",
    "Analyzing a power system for uncertain grid parameters is a relevant field for grid operation and planning.\n",
    "\n",
    "Probabilistic power flow simulations allow to efficiently analyze a grid including parameters with randomly distributed parameters. Multiple randomly distributed parameters can be convolved via Fourier transformation or FFT.\n",
    "\n",
    "Less performant but easy and flexible to use, time series simulations can produce the same results.\n",
    "This tutorial shows how random distributed parameters can be taken into account with time series simulations."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "First we start with the necessary imports."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2025-10-20T13:18:15.970555Z",
     "start_time": "2025-10-20T13:18:11.186160Z"
    },
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "from pandapower.networks.cigre_networks import create_cigre_network_mv\n",
    "from pandapower.timeseries import DFData, OutputWriter, run_timeseries\n",
    "from pandapower.control import ConstControl\n",
    "\n",
    "rng = np.random.default_rng(10)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We use the CIGRE MV grid as an example.\n",
    "We assume that both the loads at buses 4 and 9 and the generation at bus 10 are randomly distributed.\n",
    "We assume a triangular distribution for the loads and a beta distribution for the PV generation."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2025-10-20T13:18:16.710548Z",
     "start_time": "2025-10-20T13:18:15.979695Z"
    }
   },
   "outputs": [],
   "source": [
    "#  CIGRE grid\n",
    "net = create_cigre_network_mv(with_der=\"pv_wind\")\n",
    "\n",
    "# grid elements with uncertainty\n",
    "load_r4 = net.load.index[net.load.bus == 4][0]\n",
    "load_ci9 = net.load.index[net.load.bus == 9][0]\n",
    "sgen_pv10 = net.sgen.index[net.sgen.bus == 10][0]\n",
    "\n",
    "# amount of random values\n",
    "n_random_values = 1000\n",
    "\n",
    "# random values as timeseries data\n",
    "load_r4_distr = rng.triangular(0, 0.6, 1, n_random_values)\n",
    "load_ci9_distr = rng.triangular(0.3, 0.7, 1.2, n_random_values)\n",
    "sgen_pv10_distr = rng.beta(2, 2, n_random_values)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The time series data just generated looks as follows:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2025-10-20T13:18:17.499669Z",
     "start_time": "2025-10-20T13:18:16.722070Z"
    }
   },
   "outputs": [],
   "source": [
    "def hist_plot(\n",
    "        bins: int,\n",
    "        data: dict,\n",
    "    ):\n",
    "\n",
    "    colors = [\"blue\", \"orange\", \"green\", \"red\", \"purple\", \"brown\", \"pink\", \"gray\", \"olive\", \"cyan\"]\n",
    "    if len(data.keys()) > len(colors):\n",
    "        raise NotImplementedError(\n",
    "            f\"hist_plot() defines not enough colors internally. {len(colors)} colors are available \"\n",
    "            f\"but the data requires {len(data.keys())}.\")\n",
    "    fig, ax = plt.subplots()\n",
    "    for i, (key, values) in enumerate(data.items()):\n",
    "        _ = ax.hist(values, bins=bins, edgecolor=colors[i], histtype=\"step\", label=key)\n",
    "        _ = ax.hist(values  , bins=bins, color=colors[i], alpha=0.4)\n",
    "    ax.set_ylabel(\"Frequency\")\n",
    "    plt.legend()\n",
    "    plt.tight_layout()\n",
    "    return fig, ax\n",
    "\n",
    "\n",
    "fig, ax = hist_plot(50, {\"Load 4\": load_r4_distr, \"Load 9\": load_ci9_distr, \"PV 10\": sgen_pv10_distr})"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we prepare the calculation of power flows with this randomly generated input data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2025-10-20T13:18:17.547044Z",
     "start_time": "2025-10-20T13:18:17.513045Z"
    }
   },
   "outputs": [],
   "source": [
    "# construct time series DataFrame\n",
    "power_df = pd.DataFrame({\n",
    "    \"p_load_r4\"  : net.load.at[load_r4  , \"p_mw\"]   * load_r4_distr,\n",
    "    \"q_load_r4\"  : net.load.at[load_r4  , \"q_mvar\"] * load_r4_distr,\n",
    "    \"p_load_ci9\" : net.load.at[load_ci9 , \"p_mw\"]   * load_ci9_distr,\n",
    "    \"q_load_ci9\" : net.load.at[load_ci9 , \"q_mvar\"] * load_ci9_distr,\n",
    "    \"p_sgen_pv10\": net.sgen.at[sgen_pv10, \"p_mw\"]   * sgen_pv10_distr,\n",
    "    \"q_sgen_pv10\": net.sgen.at[sgen_pv10, \"q_mvar\"] * sgen_pv10_distr,\n",
    "})\n",
    "ts_data = DFData(power_df)\n",
    "\n",
    "# create constant controllers for variable parameters\n",
    "for et, var, eidx, pn in zip(\n",
    "    [\"load\", \"load\", \"sgen\", \"sgen\"],\n",
    "    [\"p_mw\", \"q_mvar\"]*2,\n",
    "    [[load_r4, load_ci9]]*2+[sgen_pv10]*2,\n",
    "    [[\"p_load_r4\", \"p_load_ci9\"], [\"q_load_r4\", \"q_load_ci9\"], \"p_sgen_pv10\", \"q_sgen_pv10\"],\n",
    "    ):\n",
    "        _ = ConstControl(net, element=et, variable=var, element_index=eidx,\n",
    "                                    profile_name=pn, data_source=ts_data)\n",
    "\n",
    "# define output writer and desired variables to be saved\n",
    "ow = OutputWriter(net)\n",
    "ow.log_variable(\"res_bus\", \"vm_pu\")\n",
    "ow.log_variable(\"res_line\", \"loading_percent\")\n",
    "ow.log_variable(\"res_trafo\", \"loading_percent\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And finally perform the power flow calculations."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2025-10-20T13:18:28.020859Z",
     "start_time": "2025-10-20T13:18:17.559889Z"
    }
   },
   "outputs": [],
   "source": [
    "# run time series\n",
    "run_timeseries(net, time_steps=range(n_random_values))\n",
    "results_dict = ow.output"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Die Netzzustände sind von den Eingangsdaten abhängig. Die resultierenden Verteilungen werden hier für die Spannungen, Transformatoren- und Leitungsauslastungen beispielhaft dargestellt:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2025-10-20T13:18:29.028299Z",
     "start_time": "2025-10-20T13:18:28.026863Z"
    }
   },
   "outputs": [],
   "source": [
    "x_label = \"Loading in percent\"\n",
    "# voltages at buses 4, 7 and 9\n",
    "bus_vm_pu = results_dict[\"res_bus.vm_pu\"]\n",
    "fig, ax = hist_plot(50, {\n",
    "    \"Bus 4\" : bus_vm_pu[4],\n",
    "    \"Bus 9\" : bus_vm_pu[9],\n",
    "    \"Bus 10\": bus_vm_pu[10]})\n",
    "ax.set_xlabel(\"Vm in pu\")\n",
    "plt.tight_layout()\n",
    "\n",
    "# loadings of lines 2, 8, 6 and 10\n",
    "line_loads = results_dict[\"res_line.loading_percent\"]\n",
    "fig, ax = hist_plot(50, {\n",
    "    \"Line 2\" : line_loads[2],\n",
    "    \"Line 8\" : line_loads[8],\n",
    "    \"Line 6\" : line_loads[6],\n",
    "    \"Line 10\": line_loads[10]})\n",
    "ax.set_xlabel(x_label)\n",
    "plt.tight_layout()\n",
    "\n",
    "# loadings of both transformers\n",
    "fig, ax = hist_plot(50, {\n",
    "    \"Trafo 0\" : results_dict[\"res_trafo.loading_percent\"][0],\n",
    "    \"Trafo 1\" : results_dict[\"res_trafo.loading_percent\"][1]})\n",
    "ax.set_xlabel(x_label)\n",
    "plt.tight_layout()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As we can see, Line 10 and Trafo 1 have no distribution since they are not effected by randomly distributed parameters located in the other feeder of the grid. Line 8 has a very narrow distribution, which can be seen more clearly in the following box plot."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2025-10-20T13:18:29.138062Z",
     "start_time": "2025-10-20T13:18:29.038201Z"
    }
   },
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots()\n",
    "results_dict[\"res_line.loading_percent\"][8].plot(kind=\"box\", vert=False, color=\"orange\", ax=ax)\n",
    "ax.set_ylabel(\"Line\")\n",
    "ax.set_xlabel(x_label)\n",
    "plt.tight_layout()"
   ]
  }
 ],
 "metadata": {
  "anaconda-cloud": {},
  "kernelspec": {
   "display_name": "Python 3.10.11 ('base')",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3"
  },
  "vscode": {
   "interpreter": {
    "hash": "19d1d53a962d236aa061289c2ac16dc8e6d9648c89fe79f459ae9a3493bc67b4"
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
